Geophysical and Astrophysical Fluid Dynamics 
Vol. 00, No. 00, October 2010, 1-18 



Effects of a- effect fluctuations on simple nonlinear dynamo models 

K. J. RICHARDSONf * and M. R. E. PROCTORf 

t DAMTP, Centre for Mathematical Sciences, University of Cambridge, Wilberforce Road, 

Cambridge, CB3 0WA 

(Accepted 3rd September 2010) 

We investigate the interaction of a fluctuating a-effect with large-scale shear in a simple nonlinear 1- 
dimensional dynamo wave model. We firstly extend the calculations of Proctor [Effects of fluctuation on 
aQ dynamo models. MNRAS 2007, 41, L39-L42] to include spatial variation of the fluctuations, and find 
that there can be a mechanism for magnetic field generation, even when the mean a is zero, provided the 
spatiotemporal spectrum of the fluctuations has an appropriate form. We investigate mean-field dynamo 
action when the new term arising from the fluctuations is non-zero, and present results concerning the 
stability and frequency of the solutions and parity selection in the nonlinear regime. The relation between 
the asymptotic theory and explicit simulation of a traditional mean-field model with a fluctuating function 
for the a-effect term is discussed. 
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1. Introduction 

It is widely accepted that the large-scale magnetic field generation in the Sun is governed by 
a dynamo process. Mathematical models of this dynamo process standardly decompose the 
magnetic field into poloidal and toroidal parts, where the aim is to complete the dynamo 
loop of generating toroidal field from poloidal field and vice versa. The presence of differential 
rotation within the Sun's convection zone is believed to be responsible for turning poloidal 
field into toroidal field, this procedure is known as the fi-effect. The mechanism by which 
poloidal field is generated from toroidal field is still open to debate, and a number of different 
models have been investigated. One such model is that of 'mean-field electrodynamics', where 
the magnetic field is supposed to exist on two very different scales, and large-scale magnetic 
fields are generated through a mean emf induced by the averaged properties of small scale 
helical motions (the 'a-effect') (see Moffatt (1978) and Krause and Radler (1980) for details). 



The a-effect is associated with an average over the small scales of the cross-product of the 
small scale velocity and magnetic fields. Since it is due to an average, it is generally treated 
as a non-zero variable varying only on long spatial and temporal scales. This assumption 
would be acceptable provided that any high-frequency fluctuations of a around this mean 
value are not too large. A turbulent, helical flow, which lacks reflectional symmetry about the 
equator, might be thought to provide perfect conditions for the a-effect to work, enhancing 
the growth of large-scale poloidal magnetic field and leading to a large-scale dynamo. A 
study by Cattaneo and Hughes (2006), however, found a contradictory result. They studied 
rotating Boussinesq convection where the flow was turbulent and helical, and found no 
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large-scale dynamo. They calculated the emf and plotted a as a function of time but found 
it to be wildly fluctuating with a very small long-time average. 

Previous work by Silant'ev (2000) and Proctor (2007) has shown that when such fluctua- 
tions in the a-effect interact with a large-scale shear, unlike Cattaneo and Hughes (2006), 
then large-scale dynamo action can occur even when the mean a is zero. Both of these 
investigations used an a-Q dynamo model; Silant'ev (2000) considered spatial variation of 
the fluctuations in quite a complicated calculation, whereas Proctor (2007) considered a 
temporally fluctuating a-effect in a much simpler setup. The Proctor (2007) calculation used 
a one-dimensional dynamo wave model and the a-effect was split into mean and fluctuating 
parts, where the fluctuations were large. A suitable average was taken, which resulted in a 
new term involving the fluctuating part of a in the mean field evolution equations. It was 
then found that this new term can lead to a large-scale dynamo mechanism, even when the 
mean a term is zero. This paper also found that when this new term was included in a simple 
mean field model of the solar cycle, the effect was a lengthening of the cycle period as the 
amplitude of the new term increased; at a sufficiently large value, the solutions became steady. 

Hughes and Proctor (2009) revisited the calculation of Cattaneo and Hughes (2006) with 
the addition of large-scale shear, and found vigorous large-scale dynamo action. They con- 
cluded that the effect of introducing large-scale shear could be due to one of two mechanisms: 
the anisotropy of the shear could lead to a shear-current effect as discussed in Rogachevskii 
and Kleeorin (2008), or the shear could interact with a temporally fluctuating a-effect as in 
Proctor (2007) both of which can be associated with off diagonal components of the turbulent 
diffusivity tensor, and lead in simple geometries to the same addition term in the mean 
field equations. Either of these effects could also explain the amplification of the large-scale 
magnetic field found in previous numerical studies by Yousef et al. (2008) and Kapyla et al. 
(2008), where 3D simulations of turbulent convection and shear were investigated. These 
calculations support the necessity of large-scale shear, and correspond to the case in Proctor 
(2007) where a large-scale dynamo mechanism was found with no mean a-effect. 

The purpose of this paper is to extend the investigation of a fluctuating a-effect in Proctor 
(2007) by looking at the linear and nonlinear consequences of the new term in a finite 
geometry. Firstly, in section 2, it is shown that the model of Proctor (2007) can be extended 
to include more realistic spatial and temporal dependent fluctuations. Detailed studies of 
the dynamics of the model in the linear and nonlinear regimes are conducted in section 3, 
which allow us to understand precisely how the steady and oscillatory modes interact, and 
the effects of parity on the system. Preliminary results of this work are discussed briefly 
in Proctor et al. (2009). Finally, in section 4, the fluctuating a-effect theory is verified by 
relating it to a simple one-dimensional model with a fluctuating function for the a-effect term. 



2. Formulation of the model 

A description of the mechanism in spherical polar coordinates can be found in Proctor (2007). 
In this paper we use a simplified one-dimensional cartesian version of the model. We begin 
with a one-dimensional Parker dynamo wave model as adapted by Proctor and Spiegel (1991). 
In this cartesian set-up, the magnetic field is assumed axisymmetric such that B = B(x, t)e y + 
V x [A(x, t)e y ], where x, y represent the North-South and azimuthal directions respectively, 
and A and B are the respective poloidal and toroidal magnetic fields. The governing equations 
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take the form of an a£l dynamo and can be written 

A t = aB + rj (A xx - I 2 A) , (la) 
B t = n'A x + r] (B xx - l 2 B) , (lb) 

where r] is the magnetic diffusivity, £ is an inverse lengthscale, and subscripts t, x denote 
differentiation with respect to time and space respectively. In equation lb, 17' is the differ- 
ential rotation (large-scale shear) which generates toroidal field from poloidal field, and the 
term involving a in equation la is the traditional a-effect term, which closes the dynamo 
loop by generating poloidal field from toroidal field by the action of small-scale helical motions. 

Spatial and temporal fluctuations are introduced in the a-effect term by writing a as a sum 
of its mean and fluctuating parts, where the fluctuating temporal and spatial scales (r and £ 
respectively) are longer than that of the original averaging process, but shorter than that of 
mean-field evolution: 

a = a + 5~ 2 ai (t,£) , (2) 

where subscripts and 1 represent the mean and fluctuating parts respectively, and depen- 
dence on x and t has been suppressed. To ensure the fluctuations are large, S 2 is taken to 
be small, and corresponds to e in the Proctor (2007) case where spatial fluctuations are not 
included. 

A method of multiple scales is applied such that 

d t ^d t + 5- 2 d T , 

d x -> d x + <5 _1 <%, 
which induces the following corrections to A and B 

A = A + A 1 (t,0 + SA 2 (t,0 + --- , (3a) 
B = B + 5B 1 (t,O + S 2 B 2 (t,C) + -- - . (3b) 

These scalings are applied to equations (la, lb), and an average is taken over the interme- 
diate scales, determined by (•} such that = (Ai) = {A 2 ) = {B{) = (B 2 ) = 0, to obtain a 
set of equations describing the evolution of the mean fields Aq and Bq: 

A 0t = a B + 5- 1 (aiBi) + ( ai B 2 ) + rj (A 0xx - £ 2 A ) , (4a) 
B 0t = n'A 0x + 77 {B 0xx - £ 2 B ) , (4b) 

where 5~ l {a.\B\) and {a\B 2 ) are new terms that arise purely from the fluctuations. The term 
<5 -1 {a\B\) appears to violate the scaling, but vanishes under plausible symmetry assumptions. 
The mean field equations in this case are given by 

A 0t = a B + ( ai B 2 ) + n (A 0xx - £ 2 A ) , (5a) 
B 0t = n'A 0x + V {B 0xx - £ 2 B ) , (5b) 

where (a\B 2 ) is the term we wish to calculate to compare with the corresponding term in 
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Proctor (2007), i.e. the term proportional to Bq x . 

In order to calculate [a\B 2 ) we subtract (5a, 5b) from the full equations (la, lb) (before 
the averaging) and take leading order terms to find equations that describe the evolution of 
the fluctuating fields A\ and By. 

A 1t = aiB + V A i^, (6a) 
B lT = n'A^ + r?Bi^, (6b) 

and we take the next order terms to find equations that describe the evolution of the fluctu- 
ating fields A2 and B 2 : 

A 2r = a x B x + 2r}A lx £ + r}A 2 ^, (7a) 
B 2t = Q'A lx + n'A 2 t: + 2r)B Xxi . + V B 2 ^. (7b) 

Equations (6a, 6b, 7a, 7b) are solved to find an expression for B 2 , which enables us to 
calculate {a\B 2 ). This average is of the form 

( ai B 2 ) = a f B - Gn'B 0x , (8) 

where a/ is a contribution to the traditional a-effect term, and 

where &\ = a\(k,uj) is the Fourier transform in space and time of a.\. There are two 
contributions to a/; one from the expansion, which vanishes under plausible symmetry 
assumptions, and the other from finite boundary effects, which are uncalculable due to the 
spatial inhomogeneity of the background field. We therefore regard a/ as a renormalisation 
of the a-effect, and concentrate on the effect of the term oc Bq x . Note that the notation 
is different from Proctor (2007), where the expression corresponding to equation (9) (with 
k = 0) is positive definite and is denoted by G 2 . The new feature of equation (9) is that when 
spatial fluctuations are included the expression may take either sign in contrast to Proctor 
(2007). More interestingly, we see that G is positive at both small and large 77, and there 
exists a range for which G is negative; (2 + ^ll/?>y l < w 2 /r? 2 fc 4 < (2 - ^fllpi)' 1 . Note 
that the corresponding formula given in Proctor et al. (2009) is not correct, although the 
nonlinear calculations are not effected by this. 

By substituting the expression for (u\B 2 ) into equation (5a), and dropping the zero sub- 
scripts on A and B, we have the final version of the model: 

A t = a B - GQ.'B x + v {A xx - £ 2 A) , (10a) 

B t = n'A x + n (b xx - e 2 B) , (iob) 

where we have subsumed a/ into the definition of ao- It is interesting to investigate 
how this new term G effects dynamo action, and whether it alone can lead to large- 
scale dynamo action with no mean a. As indicated above, other possible effects can lead 
to a term of this nature, such as the shear-current effect, see Rogachevskii and Kleeorin (2008). 
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Figure 1. Marginal stability boundary, found by plotting lines of constant m from m = 0.36 to m = 1.2 for Q against 
D 2 according to equation (11). Dynamo action is possible above the lower envelope of lines, and the region below the 
thick dashed line at Q = is the additional region of possible dynamo action found by considering spatial dependence 
in the fluctuations of a. 

The model is solved analytically by seeking solutions A, B oc exp [i (kx + tot)] to obtain the 
following dimensionless dispersion relation 
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where k is scaled such that k = £m, m constant, Q = GQ' 2 / (rj 2 £ 2 ) is the dimensionless form 
of G and D = Q'ao/(r] 2 £ 5 ) is the dynamo number. The corresponding dispersion relation 
found in Proctor (2007) has the same form but with different notation, since we now use Q 
rather than Q 2 to indicate that this term may now take either sign. 

The marginal stability boundary is plotted in figure 1 by plotting lines of constant m 
according to equation (11). Dynamo action is possible above the lower envelope of lines. This 
figure has the same form as in Proctor (2007) when Q is positive, but due to the spatial 
fluctuations, Q may be negative and the diagram can be extended as shown. It can be seen 
that dynamo action in the negative Q regime is inhibited by the new effect and larger dynamo 
numbers D are required for instability Recall that there is a range for which Q is negative, 
given above, while Q is positive when rj is both negligible and when it dominates. Figure 1 
tells us that we can achieve dynamo action when the dynamo number is zero, and therefore 
when the mean a is zero. In the case when D = 0, we find that oj = 0, and hence dynamo 
waves do not travel, so that solutions are steady in this region. To understand the stability 
of the system, and discover the range over which steady solutions can be found in a finite 
geometry, we introduce a nonlinear quenching term and solve the equations numerically. 



3. Numerical Model of the Solar Cycle 

To enable us to solve the model numerically we use the same method as Proctor (2007) to 
simulate a simple model of the solar cycle, and investigate further the stability and parity of the 
system. We introduce a nonlinear quenching so that the a-effect takes the form a(x, t)/(l+B 2 ). 
Equations (10a, 10b) are rescaled so that O' = rj = £ = 1. We also choose an. = — dsm(2irx/l) 
to simulate antisymmetric dynamo waves that travel towards the equator, where d is a positive 
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constant, and we parametrise the new term due to the fluctuations with the positive constant 
r. We are then able to vary the parameters d and r to investigate the stability in d — r phase 
space both in the linear and nonlinear regimes. The model equations then become 

—d sin (2-irx/l) B — rB x , , . , 

A t = -TTlk + Axx - A > (12a) 

B t = A x + B xx - B, (12b) 
which are to be solved between < x < I, with A = B = at x = 0, 1. 



3.1. Linear Model 

The linear stability boundary can be found by setting the denominator of the terms in B on 
the right-hand-side of (12a) to unity, leading to the system 

A t = -dsin (2irx/l) B - rB x + A xx - A, (13a) 
B t = A x + B xx - B, (13b) 

and then performing a linear stability analysis on these equations. This was done by 
expanding A and B as Fourier sine series in x and reducing the system to a matrix problem. 
This matrix was then decomposed into two individual matrices; one for dipole solutions 
(entries corresponding to A even and B odd), and one for quadrupole solutions (entries 
corresponding to A odd and B even.) The complex growth rate for each mode was then 
determined as an eigenvalue of the appropriate matrix, and the parameters d and r were 
then varied to find the stability boundaries for each parity, and for each mode of instability. 

The marginal curves for the first mode (one zero eigenvalue) and second mode (one positive 
and one zero eigenvalue) are calculated and plotted for each parity, as functions of d and r for 
varying values of Z, in figure 2. The marginal curves for the two modes (corresponding initially 
to a zero eigenvalue) join at a codimension-2 point (double zero eigenvalue), and for larger 
values of d the eigenvalues on the marginal curve are purely imaginary. The critical values of 
r at d = can be calculated exactly for this simple problem, giving r = r n = 
irrespective of parity, where n = 1 and n = 2 correspond to the first and second modes of 
stability respectively. We deduce that as I increases, r\^2 — > 4, so that the gap between the 
two modes decreases, as shown in figure 2. 

The bolder lines in figure 2 represent the marginal solution boundary; the least values of r 
for each d for marginal solutions. There are changes in parity on this boundary at the points 
where the dipole and quadrupole branches cross, however for all and large enough d, dipole 
modes are preferred at onset. The more complex parity changes are enlarged in the figures. 

All calculations have been performed for positive values of d as these are the physically 
most relevant. It is however interesting to note that due to an adjointness property of 
the governing equations the stability boundaries for negative d are the mirror images 
of those in figure 2 with dipole and quadrupole modes interchanged (Proctor 1977). It 
follows that in figure 2 the initial gradients of the marginal curves at small d for the 
two parities are equal and opposite for all values of Z. To calculate the gradient of the 
stability boundary for the first mode, r±, at d = 0, we expand A, B and r in powers of d 
in equations (13a, 13b). This results in the formula r\ = + r±^d, where is the 
value of r\ when d = 0, and r^ 1 ) is the gradient of the stability boundary at d = 0. The 
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d d 

(c) Z = 10 (d) I = 100 

Figure 2. Linear stability boundaries for dipolc (solid line) and quadrupolc (dashed line) parities, for varying values of 
I. First and second modes of stability are plotted. For sufficiently small d the eigenvalues are real. For each parity the 
curves for the two modes join at a codimcnsion-2 point, after which the eigenvalues are complex. For these two lines, the 
boundary with the lower r value at each d is shown in bold, so that the bold line gives the marginal stability boundary. 
The more complex parity changes in figures (b) and (c), and the codimension-2 point at small d in figure (d), are enlarged 
in their respective figures. 



gradient r^ 1 ) is plotted as a function of I in figure 3 for both parities. The preferred mode 
at onset corresponds to the smallest value of r\, hence quadrupole modes are preferred 
at smaller values of I, and for larger values of I the initial parity of the marginal curve oscillates. 

The frequency of the modes on the marginal solution boundary is zero before the 
codimension-2 point and increases monotonically afterwards. This corresponds to steady 
solutions initially and oscillatory solutions after the codimension-2 point. In figure 4, the 
marginal solution boundary from figure 2 is plotted (solid line), along with the corresponding 
frequency (dashed line). Vertical lines match up changes in the frequency with changes in 
parity of the marginal solution boundary. 
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Figure 3. The gradient of the marginal curve for the first mode of stability, at d = 0, plotted as a function 

of I. Dipole (solid line) and quadrupole (dashed line) parities have equal and opposite gradients. Quadrupole parity is 
preferred at small I, and the parity of the marginal curve oscillates for larger I. 




(a) 1 = 4 (b) 1 = 7 (c) I = 10 

Figure 4. The least value of r for unstable solutions for each value of d is plotted (solid line) for varying values of I. 
These lines are the bold lines from figure 2. Also plotted is the corresponding frequency (dashed line) for each point on 
the line. Changes in the parity of the marginal solution boundary align with discontinuities in the frequency. Points on 
the marginal solution boundary before the codimension-2 point have zero frequency. 



3.2. Nonlinear Model 

To investigate the behaviour of the system in the nonlinear regime, we consider the full 
nonlinear system given by (12a, 12b) which include the simplest form of a quenching for 
both the mean and fluctuating parts of the a-effect. It transpires later that the direct 
application of the fluctuating theory to this nonlinear system does in fact lead to a different 
nonlinearity in the r term, see (18a, 18b), however, solutions with this different nonlinearity 
are qualitatively much the same. 



These equations are solved using a second order finite difference scheme, in a box of size 
< x < I, where x = represents the North pole, and x = I represents the South pole. This 
allows us to define the two parities such that dipole parity corresponds to B odd and A even 
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Figure 5. (a) Period of oscillations plotted against r (the amplitude of the fluctuating a), for varying values of d (the 
amplitude of the mean a). The values of d are d = 5 (solid line), d = 10 (dashed line), and d = 15 (dotted line), (b) 
Regions of different stability in d-r phase space. 



about the equator (x = 1/2), and quadrupole parity corresponds to B even and A odd. A and 
B are set to be zero at the boundaries, and d and r are varied to investigate stability and parity 
changes in d—r parameter space. In all of the nonlinear investigations a value of I = 10 is used. 

Proctor (2007) found that for a fixed value of d = 30, and varying values of r, both 
oscillating and steady solutions can be found. In fact, the period of the oscillating solutions 
increases with increasing values of r, and steady solutions are achieved only when r is 
sufficiently large. Further investigation reveals that varying the value of d, the critical value 
of r to achieve steady solutions, r c , also changes. Furthermore, the larger the value of d, the 
larger the value of r c . This is shown in figure 5(a) where the period of oscillations is plotted 
as a function of r for three different values of d. In this figure, a steady solution is defined as 
one with infinite period. The increase of r c with d is clear. 

The onset of steady solutions in d — r space is plotted in figure 5(b). The marginal stability 
boundary for I = 10 from the linear model has also been included in the figure and the regions 
of different types of solution are labelled. As r increases through the oscillating region, the 
period of the solutions increases until it becomes infinite and the solution reaches a steady 
state. 

Mixed initial conditions were used in figure 5(b) so that neither A nor B are even or odd 
about the equator. This allows solutions to settle into their preferred parity after integrating 
for long times. A parameter study reveals that there are many different regions of stable 
parity throughout parameter space. These are shown in figure 6(a), where the boxes are 
enlarged in figures 6(b) and 6(c) to show the more detailed structure. There are large regions 
of mixed parity, where the solution does not settle into either dipole or quadrupole parity 
after integrating for long times, between regions of dipole and quadrupole parity in the 
steady regime. However, the only region of mixed solutions in the oscillating region is that in 
figure 6(c) separating the initial parity from the stability boundaries. We see that this region 
of mixed solutions meets the stability boundary at exactly the same point as where the two 
separate parity stability boundaries cross (seen in figure 2). 
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(b) (c) 

Figure 6. (a) Regions of stable parity throughout phase space. The two boxes in figure (a) arc enlarged in figures (b) 
and (c) below. 

It is possible to fix the solutions to be either dipole or quadrupole by integrating equations 
(12a, 12b) over half of the domain. It is then necessary to fix the boundary conditions so that 

A = at x = 0, 

B x = at x = l/2, 

for a dipole solution, and 

A x = at x = 0, 

B = at x = 1/2, 



for a quadrupole solution. 
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Figure 7. Regions of different stability for fixed dipole (solid) and quadrupole (dashed) parities. The onset of steady 
dynamo action is plotted for each parity and the stability boundaries from figure 2(c) are added. The region where the 
steady boundary meets the stability boundary at the codeimension-2 point for each parity is enlarged. 



Figure 7 shows the onset of steady solutions for the two separate parities, obtained by 
using these boundary conditions. The critical value of r to achieve steady solutions is lower 
for quadrupole parity solutions than it is for dipole parity solutions. In this figure, the 
stability boundaries for both dipole and quadrupole parities from figure 2(c) are added. For 
each parity, the line plotting the onset of steady solutions meets the stability boundary at 
the codimension-2 point, which is enlarged in the figure. This reinforces the fact that steady 
solutions are found to the left of this point, and oscillatory solutions are to the right. 



4. Testing the asymptotic model with a rapidly varying a function 

The theory leading to the model we have studied can be tested against a simple one- 
dimensional dynamo model (similar to (la, lb) adapted to be solved numerically in the 
same way as in section 3) where a is now a fluctuating function of time with mean value d 
and fluctuating part a. By increasing the magnitude of a, we expect to see the cycle period 
increase, and if the theory is correct there should be a predictable relation between a, r and 
the increase in period. 

The model used to test the theory is 

A t ^- i ^lf B t S -^ B + A,,-A, ,14a) 
B t = A x + B xx - B, (14b) 

where d is the same mean a term as in the fluctuating theory, 5 is a large parameter to 
ensure a is large (the same scaling as is used in section 2), and saturation is included so 
that these equations can be solved numerically, and to ensure the magnetic fields are quenched. 
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The form of a is given by 



a = ( a i s'm(iat), 



(15) 



where a is the frequency, which is large such that a = 5~ 2 a where a is on order 1 quantity, 
( is a parameter used to vary the amplitude of the function, and CKj is a sample of randomly 
generated numbers from a normal distribution with unit standard deviation. A time series 
for a is plotted in figure 8. 

To be able to understand the relationship between a and r in our fluctuating theory, 
we recall from Proctor (2007), where a time-dependent fluctuating a-effect was used, that 
{ot\Bi) = —G 2 Q'B X} and recalculate (a\Bi) using equation (15) for oi\. Recalling that when 
transforming to the numerical model, ffl = r\ = 1 and G 2 — > r, we find that r and £ are related 
as follows 

i=i 

This expression predicts the corresponding values of r and £ that should give solutions with 
the same period. Increasing the value of r or ( according to this expression should result in 
the same increase in period. 



4.1. New form of nonlinear fluctuating ol model 

In order to provide a quantitative comparison between the models in the nonlinear regime it 
is necessary to derive a new more accurate set of equations for the fluctuating theory. This 
is because the new model, (14a, 14b), possesses saturation, whereas saturation is only added 
to the fluctuating theory when adapting it to the numerical model. To be able to compare 
these two models we must rederive the fluctuating theory equations from first principles, with 
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Tabic 1. A comparison of the period for the asymptotic model for 
selected values of r with the period of the test model with corre- 
sponding £ (calculated according to equation (16)). 



r 


Period for asymptotic theory 


Period for test model 





1.53 


1.53 


0.50 


1.63 


1.63 


2 


1.86 


1.88 


2.02 


1.88 


1.85 


4 


2.18 


2.17 


4.54 


2.28 


2.30 


6 


2.52 


2.42 


8 


2.82 


2.73 


8.07 


2.84 


2.73 


10 


3.13 


3.00 


12 


3.40 


3.15 



a quenched a-effect. The governing equations are then (compare with (12a, 12b)): 

A t = (a + 5' 2 ai )f{B) + rj{A xx - I 2 A), (17a) 
B t = n'A x + r l (B xx -e 2 B), (17b) 

where ao is the mean part of a, ot\ is the fluctuating part of a and is a function of the small 
timescale r, and f(B) = B/(l + B 2 ) for our model. 

Following the same method as in sections 2 and 3, a new nonlinear model can be derived 
(see the Appendix), given by 

dsm(27rx/l)B r(l-B 2 ) 2 B x , , x 

A t = v ' * y - + A xx - A, 18a 

1 + B 2 (1 + B 2 ) 4 y ' 

B t = A x + B xx - B. (18b) 
The solutions of this model are qualitatively similar to the original nonlinear model. 



4.2. Validating the asymptotic theory 

The test model with a fluctuating a function, (14a, 14b), and the new nonlinear version of 
the fluctuating theory, (18a, 18b), are both solved in the same way as described in section 
3.2. We set I = 50, d = 25 for both models so that when both the new term in the asymptotic 
model (r) and the new fluctuating function (a) are zero, both models will have the same 
solution with the same period. We use a = 100 in the test model for rapid temporal variations 
in a, and choose 5 = 0.1 so that a = 1. 

The analytic prediction, equation (16), is tested by solving the improved asymptotic 
model (18a, 18b) with a chosen value of r, and solving the test model (14a, 14b) with the 
corresponding value of according to equation (16), and the periods of both solutions are 
measured and compared. Results are shown in Table 1, where some entries correspond to 
choosing ( first and then calculating the corresponding value of r. 

The results in Table 1 are plotted in figure 9(a). The solid line represents the periods from 
the improved nonlinear model, and the results from the test model are plotted as crosses. It 
is clear that when r is small, the models are very similar, however, the difference between 
the periods gets larger as the solution moves towards a steady state at larger values of r. In 
figure 9(b), for r = 10, where the difference in period between the two models is sufficiently 
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Figure 9. (a) A comparison of r against period for the improved nonlinear model (solid line) and the test model (crosses), 
using data from Table 1. The periods are very similar for small values of r, but as r increases, the difference between the 
periods for the two models also increases, (b) A comparison of the toroidal field B plotted against x for the improved 
nonlinear model (solid line) and the test model (dashed line) for r = 10 and ( = 4.45. 

large, the toroidal field, B, is plotted against x. The improved nonlinear model (solid line) 
and the test model (dashed line) are very much the same, however they do differ slightly, 
hence the difference in period of these solutions. 

A possible explanation for the difference between the models at large r is that the asymptotic 
theory doesn't allow for the fact that A and B have to vanish at the boundaries, so there 
are likely to be small errors associated with this. In an attempt to reduce this error, we have 
used a larger value of Z, Z = 50, which has minimised the boundary effects. However, the 
general conclusion is that the overall effect of increasing the amplitude of a in the test model 
is an increase in the period of the solution, and therefore the test model seems to correctly 
reproduce the results from the time-dependent asymptotic a theory at small r. 

4.3. Spatiotemporal Case 

Since one of the aims of this paper is to investigate the addition of spatial dependence in the en- 
effect fluctuations, it is of interest to extend the verification of the theory to the spatiotemporal 
case. We now take a to be 

a = C, oiij sin (iat) sin (jnx) , (19) 

i=l j=l 

where £, a^j and a are the same as in equation (15), and k = 5~ l R is the spatial fequency, 
where k is an order 1 quantity. 

Unfortunately, a quantitative comparison with the spatiotemporal asymptotic theory is 
not practical because it is not possible to calculate the contribution to the a-effect in a 
finite geometry (as mentioned in section 2). However, it is possible to solve equations (14a, 
14b) with a as given in equation (19), in the linear regime (for simplicity), and show that 
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Figure 10. f 2 plotted as a function of the growth rate of magnetic energy for the test model with a spatially and 
temporally fluctuating a-effect. 

increasing £ does indeed increase the growth rate of magnetic energy. 

The theory allows arbitrary additive choices of different modes, but it is simpler to use just 
one mode for calculations, due to numerical difficulties with resolution, therefore we take 
imax = jmax = 1- We set Z = 50, <i = 25 and choose a = 1 and R = 0.5; this corresponds 
to a positive value of r, where we would expect to see an increase in the growth rate in the 
asymptotic theory. Results are shown in figure 10, where Q 2 is plotted against the growth 
rate. Although it is not possible to make a comparison between the two models in the 
spatiotemporal case, it is clear that an increase in ( does lead to an increase in the growth 
rate of magnetic energy, as predicted by the asymptotic theory. 



5. Conclusion 

The interaction of a spatially and temporally fluctuating a-effect with large-scale shear can 
lead to a large-scale dynamo mechanism, as in Proctor (2007). Moreover, a similar term 
to the new term found in Proctor (2007) relating to the fluctuating a-effect is calculated 
and found to be positive or negative, depending on the spatiotemporal spectrum of the 
fluctuations. If the term is negative the critical value of the dynamo number is increased, and 
dynamo action is inhibited. 

Numerical investigations were made, both in the linear and nonlinear regime, by adapting 
the equations to a simple one-dimensional model of the solar cycle, and investigating 
behaviour in the two-dimensional d — r parameter space. In the linear model it was found 
that the values of r for the onset of unstable solutions are different for dipole and quadrupole 
parity systems at a fixed value of d, but these become more and more alike as I increases. The 
first two modes of stability are initially steady state, but join together at a codimension-2 
point, where there is a region of solutions with zero frequency, i.e. steady solutions, to the 
left of this point, and oscillating solutions to the right. 

In the nonlinear model, the onset of steady solutions was plotted in d — r space and it 
was found that an increase in the fluctuating a parameter, r, led to an increase in the 
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cycle period of the oscillations. Regions of stable parity were determined throughout phase 
space and were found to be quite complicated considering the simplicity of the model. The 
boundary conditions were adjusted to allow the determination of the onset of steady and 
oscillatory solutions for the dipole and quadrupole cases separately, and combining these 
with the corresponding stability boundaries from the linear model reveals that the curve 
separating steady solutions from oscillatory solutions in the nonlinear regime tends to the 
codimension-2 point where the two types of linear stability join. 

Finally, the asymptotic theory is verified against a one-dimensional dynamo wave model 
with a fluctuating function for the a-effect. Firstly, the time-dependent case is verified in 
the nonlinear regime, where a new more accurate version of the asymptotic theory from 
Proctor (2007) is derived and used to compare with a temporally fluctuating function for 
a. The relationship between the time-dependent function and the previous r term is found 
and tested, and we see an increase in the cycle period of the oscillations in the predicted 
way according to the relationship derived. Secondly, an attempt to verify the spatiotemporal 
case is made using a space and time-dependent a function. It is not possible to make a 
quantitative comparison with the asymptotic theory in this case, due to contributions to the 
a-effect in the analytic theory which cannot be calculated in a finite domain. Instead, the 
test model is solved, in the linear regime, and it is shown that increasing the value of the 
spatiotemporal a function does indeed increase the growth rate of magnetic energy. 



Appendix A. 

Derivation of improved nonlinear model 

The governing equations are 

A t = af(B)+ v (A xx -£ 2 A), 
B t = n'A x + r](B xx -£ 2 B). 

We follow the derivation of the fluctuating a model used in Proctor (2007), since this is the 
temporal case, with a more general function of B. We split a, A and B into their mean and 
fluctuating parts and apply a method of multiple scales such that 

a = a + e _1 ai(r), 
A = A + A 1 (t), 
S = 5 + eB 1 (r), 
f(B) = f(B ) + eB 1 f'(B ) + -- - , 

where subscripts and 1 denote the mean and fluctuating parts respectively, r is the in- 
termediate timescale and e is taken to be small. A suitable average is taken such that 
(ai) = (Al) = (Bi) = and (f(B)) = f(B ) + 0(e 2 ), and we obtain the following lead- 
ing order mean field equations 



A 0t = a f(B ) + (aiBi) f'(B ) + v (A 0xx - £ 2 A ), 
B 0t = n'A 0x + ri(B 0xx -£ 2 B ). 



(Ala) 
(Alb) 



Effects of fluctuations on nonlinear dynamo models 



17 



The new term {a\B\) is calculated in the usual way and found to be 

(aiSi) = -G 2 f'(B )n'B 0x . 

Substituting (a\B\) in Equation (Ala) and, to use the same notation as section 3.2 for 
comparison, we set ao = —dsm(2irx/l), G 2 = r,ft' = ri = £= l and drop the zero subscripts. 
Equations (Ala, Alb) then become 

A t = -dsm (2ttx/1) f{B) - r(f'(B )) 2 B x + A xx - A, (A2a) 
B t = A x + B xx - B. (A2b) 

In section 4.1, f(B) = B/(l + B 2 ), so we find 

f'(B)= l ~ B2 



(1 + S 2 ) 2 ' 
and equations (A2a, A2b) become 

dsva.{?irx/l)B r{l-B 2 ) 2 B x 

Bt = A x + B xx - B. (A3b) 

This derivation only concerns a temporally fluctuating a-effect. It would be impracticable 
to calculate the corresponding improved nonlinear model in the spatiotemporal case due to 
the multiple levels of expansion. 
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